Drag Reduction by Polymer Additives in Decaying Turbulence 
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We present results from a systematic numerical study of decaying turbulence in a dilute polymer 
solution by using a shell-model version of the FENE-P equations. Our study leads to an appealing 
definition of drag reduction for the case of decaying turbulence. We exhibit several new results, such 
as the potential-energy spectrum of the polymer, hitherto unobserved features in the temporal evo- 
lution of the kinetic-energy spectrum, and characterize intermittency in such systems. We compare 
our results with the GOY shell model for fluid turbulence. 
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The phenomenon of drag reduction by polymer 
additives^, whereby dilute solutions of linear, flexible, 
high-molecular-weight polymers exhibit frictional resis- 
tance to flow much lower than that of the pure sol- 
vent, has almost exclusively been studied within the con- 
text of statistically steady turbulent flows since the pi- 
oneering work of TomsQ- By contrast, there is an ex- 
treme scarcity of results concerning the effects of poly- 
mer additives on decaying turbulence 3]. Experimen- 
tal studies of decaying, homogeneous turbulence behind 
a grid indicate, for such dilute polymer solutions, a 
turbulent energy spectrum similar to that found with- 
out polymers jj, |5j- However, flow visualization via 
die-injection tracers^ and particle image velocimetryj^ 
show an inhibition of small-scale structures in the pres- 
ence of polymer additives. To the best of our knowl- 
edge decaying turbulence in such polymer solutions has 
not been studied numerically. We initiate such a study 
here by using a shell model that is well suited to ex- 
amining the effects of polymer additives in turbulent 
flows that are homogeneous and in which bounding walls 
have no direct role. We obtain several interesting results 
including a natural definition of the percentage drag- 
reduction DR, which has been lacking for the case of de- 
caying turbulence. We show that the dependence of DR 
on the polymer concentration c is in qualitative accord 
with experiments as is the suppression of small-scale 
structures which we quantify by obtaining the filtered- 
wavenumber-dependence of the flatness of the velocity 
field. 

We will use a shell-model version of the FENE-P (Finitely 
Extensible Nonlinear Elastic - Peter lin) 0,0 model for di- 
lute polymer solutions that has often been used for study- 
ing viscoelastic effects since it contains the basic charac- 
teristics of molecular stretching, orientation and finite ex- 
tensibility seen in polymer molecules. A direct numerical 
simulation of the FENE-P equations is computationally 
prohibitive. This motivates the use of a shell model that 
captures the essential features of the FENE-P equations. 
Recent studies [tj have exploited a formal analogy ^(J of 
the FENE-P equations with those of magnetohydrody- 



namics (MHD) to construct such a shell model. We in- 
vestigate decaying turbulence in a dilute polymer solu- 
tion by developing a similar shell model for the FENE-P 
equations. The unforced FENE-P equations 0,0] are 
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where p is the pressure, v s the kinematic viscosity of the 
solvent, v p a 'viscosity' parameter, p s the density of the 
solvent, incompressibility is enforced via sj ■ v = 0, and 
the polymer conformation tensor is H a p = {RaRp) / 'Rq, 
with the angular brackets indicating an average over 
polymer configurations, of the dyadic product of the end- 
to-end vector R(x,t) of the polymer molecules. The max- 
imal extension of the polymer molecules is restricted by 
the condition (i? 2 ) < Rq. The contribution to the stress 
tensor because of the polymer is T a p — v p [P(x,t)lZ a p — 
5 a p]/T p , with S a /3 the Kronecker delta, r p the time con- 
stant of the FENE-P model, and P(x, t) = 1/(1 - ft 77 ) 
(with repeated indices indicating a trace). The concen- 
tration of the polymer is parametrized here by c = v p /u s . 
Our shell-model version of the unforced FENE-P equa- 
tions, obtained by generalising a shell model originally 
proposed for three-dimensional MHD 11], is 

^ = ^n.w - v s k 2 n v n + —P{b)<f> n , bb , 
dt t„ 



db n 
dt 



$n,vb ~ ®n,bv ~ —P(b)b ri 



(2) 



where P{b) = 1/(1 — X] n |6„| 2 ), v n and b n are com- 
plex, scalar variables representing the velocity and 
the (normalized) polymer end-to-end vector fields, 
respectively, with the discrete wavenumbers k n = ko2 n 
(ko sets the scale for wave- numbers), for shell in- 
dex n (n — 1...N, for AT shells), with $ n , vv = 

i(aik n V n+ lV n +2 + a 2 kn-lV n+ lV n -l + a 3 k n -2V n -lV n -2), 



2 



®n,bb = -i{aik n b n+1 b n+2 + a 2 k n - 1 b n+1 b n - 1 + 
a3k n -2b n -ib n _ 2 ), $n,vb = i{a4,k n v n+ ib n+ 2 + 
azkn-iVn-iK+i + a 6 k n -2V n -iK-2) , and § nM = 
-i(a 4 k 

nbn+l^jjA-2 ~t~ ^5^n— l^n— l^n+1 ^6^n — 2^n — l^n — 2)- 

As in Ref. [ljj we choose a\ = 1, a 2 = (13 = —1/2, 
a 4 = 1/6, 05 = 1/3, and a§ = —2/3. We solve 
Eqs. (0 numerically by using an Adams-Bashforth 
scheme 13] and double-precision arithmetic, with a step 
size 5t = 1(T 2 and N = 22 shells, with fc = 1/16 and 
v s = 10~ 3 . For all our runs (except those in FIGs. Eland 
0} we set c = 100. For numerical stability we add a nom- 
inal viscous term —v^k^bn to the shell-model equations 
for b n and set Vb/v s = 10~ 13 . With these parameter 
values, our code is stable for 1.0 < t p < 7.8, and we ob- 
serve that the corresponding percentage drag reduction 
DR (see below) lies in the range 63% < DR < 98%. For 
specificity we use t p — 2.1 for the data presented here. 
The initial velocity field is taken to be V 7 -. - fan G n 
(for n = 1,2), v° = kl/ 2 e~ k ^e ie " (for 3 < n < N) 
and the initial polymer field to be 6^ ~ fc„ 1 / 2 e l< ^ > ", 
with 6 n and </>„ independent random phases distributed 
uniformly between and 2ir. In decaying turbulence, 
it is convenient to measure time in units of the initial 
large eddy-turnover time. For our shell model this 
is r = l/{v° rms k x ) with v° rms = [(E„l<| 2 )] 1/2 , ^e 
root-mean-square value of the initial velocity (we find 
To = 5.2). We use the dimensionless time r = t/ro (t 
is the product of the number of steps and St). Our 
runs are ensemble averaged over 10 4 independent initial 
conditions with different realizations of phases. We 
define Re® = v^ ms / {k\v s ) to be the value of the initial 
Reynolds number (here Re® equals 12309). Shell- model 
energy densities are defined as EAkjj) = {\a n \ 2 /k n ), 
with a = v for the velocity fieldyjj, [l2j and a — b for 
the polymer field. Equations @ reduce to those for 
the Gledzer-Ohkitani-Yamada (GOY) shell model[I3 
when the polymer-field terms are suppressed. For our 
GOY shell model runs, we use initial parameter values 
as for the FENE-P shell model to facilitate comparisons 
between the two. 
Figure ^a) shows the time evolution of the normalized 
kinetic energy spectrum E v {k n ) / 'E v {k\) (successive 
curves separated by time intervals of 0.2r). We see a 
cascade of the energy to large wavenumbers after which 
the shape of the spectrum does not change appreciably 
but the energy decays. We observe the evolution of 
a flat portion in the spectrum that vanishes upon 
cascade completion (plot with open circles). Figure ^b) 
compares kinetic-energy spectra at cascade completion 
for our model (Eqs. J2Jl) and for the GOY shell model. 
In the inertial range, both spectra are indistinguishable 
and show a Kolmogorov-type fc -5 / 3 behavior with an 
observed slope of —1.67 ± 0.01 (with errors from least- 
square fits), a result consistent with experimentsyy,|5j of 
decaying, homogeneous turbulence behind a grid for a 
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FIG. 1: (a) Log-log plots of the temporal evolution of the nor- 
malized kinetic energy spectra E v (k„) /E v (ki) of the FENE-P 
shell model as a function of the wavenumber k n . The plot with 
open circles is calculated at cascade completion, 
(b) Log-log plots of the normalized kinetic energy spectra 
E v (k n )/ E v (ki) as a function of the wavenumber k n for the 
FENE-P and the GOY shell models at cascade completion. 
The observed slope is —1.67 ±0.01 for the range 0.25 < k n < 
64. 



dilute polymer solution. However, significant differences 
show up in the dissipation range: the spectrum for 
the FENE-P shell model falls much more slowly than 
its GOY-model counterpart indicating greatly reduced 
dissipation at large wavenumbers. Experimental^, Q 
energy spectra do not cover as large a range of spatial 
scales as we can cover in our shell-model study, and thus, 
to the best of our knowledge, these dissipation-range 
discrepancies of the energy spectra, with and without 
polymer additives, have not been noticed earlier. We 
note that our results in FIG. db) distinctly differ from 
corresponding results for statistically steady turbu- 
lence, where a tilt in the spectrum has been observed at 
low wavenumbers in the FENE-P shell model relative to 
that obtained from the GOY shell model. 

In FIG. Ela), we display the time evolution of the 
potential-energy spectrum of the polymer Eb(k n ) (with 
a temporal separation of 0.2r). Starting from an initially 
flat spectrum, we observe the appearance and subsequent 
growth of a protuberance that bulges out maximally 
on cascade completion (plot with open circles) at a 
wavenumber corresponding to the value, of order unity, 
of the ratio of the polymer time constant r p and the 
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FIG. 2: (a) Log-log plots of the temporal evolution of the 
potential energy spectrum Eb(k n ) of the polymer field as a 
function of the wavenumber k n . The plot with open circles is 
calculated at cascade completion. 

(b) Log-log plot of the ratio of the time constant r p of the 
FENE-P shell model and the turbulence time-scale r(k„) as 
a function of the wavenumber k n at cascade completion. The 
inset shows the inverse of the turbulence time-scale r(fe n ) as 
a function of the wavenumber k n for the GOY shell model at 
cascade completion. 
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FIG. 3: The normalized kinetic energy dissipation rate 
£f/£f,o °f the FENE-P shell model as a function of the di- 
mensionless time r for different values of concentration c, as 
specified in the legend. 
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FIG. 4: Plot of the percentage drag reduction DR (Eq. 
as a function of the concentration c. The inset shows a plot 
for the same quantities taken from Ref. [l5| (see the text for 
definitions). 



turbulence time scale r(fc„) = l/(k n \v(k n )\) (FIG. 
Elb)). The result is in agreement with a hypothesis (for 
statistically steady turbulence) in Ref. |14j wherein a 
polymer molecule, immersed in an eddy with a turbulent 
time-scale comparable to the polymer relaxation time, 
undergoes a 'coil-stretch' transition with an increment 
in the potential-energy spectrum at the wavenumber 
corresponding to the inverse of the eddy size. The inset 
in FIG. |2Ib) is a plot of the inverse of the turbulence 
time-scale r(k n ) as a function of the wavenumber k n for 
the GOY shell model. In both plots, within the inertial 
range, r(fc„) ~ fc„ 2 ^ 3 , a result consistent with the —5/3 
power-law in the kinetic energy spectrum. 
A log-log plot of the normalized kinetic energy dissi- 
pation rate £f/£f t o of the FENE-P shell model versus 
dimensionless time r for different values of c is shown in 
FIG. [3 (with £f/ £ffi = (£ n k*\v n \ 2 )/(£ n k 2 n \v° n | 2 >, the 
additional index indicating values calculated at initial 
times). The reduction in the peak value with respect to 
the value at initial times, with increasing concentration, 
is indicative of an enhanced value of the stored elastic 



potential energy in the polymer molecules due to their 
extension. The analogous plot for the GOY shell model 
is identical to the plot for c = 10 in FIG. We are 
therefore led to the following natural definition of the 
percentage drag reduction DR for decaying turbulence: 

DR = f £ 3^/£g,0-£f,rn/£f,0 \ x ^ ^ 

where the kinetic energy dissipation rates £ a = 
(En k lK\ 2 ) (the subscript a = f for the FENE-P shell 
model and a = g for the GOY shell model) are calcu- 
lated upon cascade completion when the dissipation rate 
is a maximum (indicated by an additional subscript m) 
and normalized by their values at initial times (indicated 
by an additional subscript 0). With the choice of ini- 
tial parameter values as specified above, £ g . m /£ g ,a equals 
234.96 ±0.01. 

In FIG. 01 we use Eq. @ to plot DR as a function of c. 
The inset figure from Ref. 1151 is a similar plot for a di- 
lute solution of Carrageenan[ljj| (a seaweed derivative) in 
a pipe-flow Reynolds number of 14000. The qualitative 
agreement with a laboratory experiment (for statistically 
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FIG. 5: Plot of the flatness Fk as a function of the filtered 
wavenumber K (see the text for definitions) at cascade com- 
pletion for the FENE-P and the GOY shell models. 



steady turbulence) supports our definition of drag reduc- 
tion for decaying turbulence. 
Laboratory experiments in both statistically steady [l^ 
and decaying|]J turbulent flows of dilute polymer 
solutions show an inhibition of small-scale structures 
and narrower probability distribution functions of ve- 
locity differences fl7j. We plot in FIG. [5] the flatness 
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(E„ l u n| 2 ), with n = K...N), K = 1...N, as a function of 
the 'filtered' wavenumber K for the FENE-P and the 
GOY shell models at cascade completion. We observe 
that, in the GOY shell model, the flatness Fk exhibits 
unbounded growth for large wavenumbers, an indication 
of strong intermittency in the dissipation scales. How- 
ever, for the FENE-P shell model, we observe that the 
flatness is greatly reduced relative to that for GOY and, 
in fact, decreases in the dissipation scales. Our results are 
consistent, therefore, with laboratory experiments which 
show a suppression of small structures that would im- 
ply reduced intermittency in the dissipation range of our 
shell model. 

Laboratory experiments [f| 0] of decaying turbulence be- 
hind grids indicate a reduced decay rate of the kinetic 
energy in a dilute polymer solution, relative to the pure 
solvent. In the initial period of decay, before the integral 
scale of turbulence becomes of the order of the size of the 
system (the minimum wavenumber, in the case of shell 
models), we observe a decay rate of —1.80 ± 0.01 for the 
FENE-P shell model and a decay rate of -2.01 ± 0.02 for 
the GOY shell model (a result consistent with Ref. (l^h 
In conclusion, then, we have presented results from a 
systematic numerical study of decaying turbulence in a 
dilute polymer solution by employing a shell-model for 
the FENE-P equations. This leads to a natural defini- 
tion of drag reduction for such a system and new results 
on the potential- and kinetic energy spectra which are in 
qualitative agreement with experimental findings. 
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